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A simple yet general method for constructing basis sets for molecular electronic structure calcu- 
lations is presented. These basis sets consist of atomic natural orbitals from a multi-configurational 
self-consistent field calculation supplemented with primitive functions, chosen such that the asymp- 
totics are appropriate for the potential of the system. Primitives are optimized for the homonuclear 
diatomic molecule to produce a balanced basis set. Two general features that facilitate this basis 
construction are demonstrated. First, weak coupling exists between the optimal exponents of prim- 
itives with different angular momenta. Second, the optimal primitive exponents for a chosen system 
depend weakly on the particular level of theory employed for optimization. The explicit case con- 
sidered here is a basis set appropriate for the Burkatzki-Filippi-Dolg pseudopotentials. Since these 
pseudopotentials are finite at nuclei and have a Coulomb tail, the recently proposed Gauss-Slater 
functions are the appropriate primitives. Double- and triple-zeta bases are developed for elements 
hydrogen through argon. These new bases offer significant gains over the corresponding Burkatzki- 
Filippi-Dolg bases at various levels of theory. Using a Gaussian expansion of the basis functions, 
these bases can be employed in any electronic structure method. Quantum Monte Carlo provides 
an added benefit: expansions are unnecessary since the integrals are evaluated numerically. 



I. INTRODUCTION 

In quantum chemistry (QC) calculations, molecular 
orbitals are traditionally expanded in a combination of 
primitive Gaussian basis functions and linear combina- 
tions of Gaussian primitives called contracted basis func- 
tions t!|. These basis sets cannot express the correct 
molecular orbital asymptotic behavior but are used in 
QC calculations to permit analytic evaluation of the two- 
electron integrals [2] . 

Analytic integral evaluation significantly limits flexibil- 
ity in basis set choice but is essential for computational 
efficiency in QC calculations. However, in practice, other 
basis function forms can be considered since an arbitrary 
function can be expanded in Gaussians. Of course, the 
fidelity of this representation is limited. An expansion in 
a finite number of Gaussians cannot reproduce the expo- 
nential decay of the wavefunction at large distances or 
the Kato cusp conditions [3J at nuclei, but it can mimic 
these features over a finite range. 

Quantum Monte Carlo (QMC) calculations [j] offer 
greater freedom in choice of basis functions because ma- 
trix elements are evaluated using Monte Carlo integra- 
tion. Consequently, the correct short- and long-distance 
asymptotics can be satisfied exactly. For systems with 
a divergent nuclear potential. Slater basis functions can 
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exactly reproduce the correct electron-nucleus cusp and 
long-range asymptotic behavior of the orbitals. For cal- 
culations on systems with a potential that is finite at 
the nucleus and has a Coulomb tail, Gauss-Slater (GS) 
primitives I5| are the appropriate choice since they intro- 
duce no cusp at the origin and reproduce the exponential 
long-range asymptotic behavior of the orbitals. 

Despite shortcomings, traditional QC basis sets have 
yielded good results. The natural orbitals (NOs) from a 
post Hartree-Fock (HF) method are aparticularly suc- 
cessful form of contracted function p-Q. The simplest 
NO construction involves diagonalizing the one-particle 
density matrix from a ground state atomic calculation 
6] . This construction is unbalanced due to obvious bias 
favoring the atom. More complicated constructions in- 
volve diagonalizing the average one-particle density ma- 
trix of several systems: atomic ground and excited states, 
ions, diatomic molecules, and atoms in an external elec- 
tric field ^-i9|]. These constructions produce excellent 
results, but they are complex. 

A simple but general method for constructing basis sets 
for molecular electronic structure calculations is proposed 
and tested here. The bases are combinations of the NOs 
obtained from diagonalizing the one-particle density ma- 
trix from an atomic multiconfigurational self-consistent 
field (MCSCF) calculation and primitive functions ap- 
propriate for the potential in the system. The primi- 
tives are optimized for the homonuclear dimer in cou- 
pled cluster calculations with single and double excita- 
tions (CCSD), with the intention of producing a balanced 
basis set. Importantly, optimal exponents for the prim- 
itive functions are shown to depend weakly on the level 



of theory used in the optimization. Additionally, results 
show that coupling is weak between primitive functions 
of different angular momenta. This enables efficient de- 
termination of optimal exponents. 

The utility of the above construction is demonstrated 
for the elements hydrogen through argon with the non- 
divergent pseudopotentials of Burkatzki, Filippi, and 
Dolg (BFD) [10]. Since these pseudopotentials are fi- 
nite at the nuclei and have a Coulomb tail, the GS 
functions are the appropriate primitives. These pseu- 
dopotentials are chosen for demonstrated accuracy in all 
cases tested and because they are accompanied by a ba- 
sis set. The BFD basis [l3| serves as a metric for test- 
ing the new basis. The benefits of our bases extend to 
all electronic structure methods tested, including CCSD, 
HF, the Becke three-parameter hybrid density functional 
(B3LYP) 11], and QMC. 

The main area of interest for the authors is QMC. Since 
QMC results depend less on basis set than traditional QC 
methods 0, only double-zeta (2z) and triple-zeta (Sz) 
bases are presented. 

This paper is organized as follows. Basis function form 
and properties are demonstrated in Sec. |lTl Results for 
calculations with the new bases are discussed in Sec. IIIII 
Concluding remarks are provided in Sec. IIVI Supple- 
mentary material is provided on EPAPS [12l ]. 



II. BASIS SET 

The number of basis functions for each angular mo- 
mentum follows the correlation consistent polarized basis 
set prescription of Dunning [1| . 2z and 3z bases appro- 
priate for the BFD pseudopotentials are generated for 
the elements hydrogen through argon. Since the BFD 
pseudopotential removes no core for hydrogen and he- 
lium, the 2z basis for these elements consists of two S 
functions and one P function, while the 3z basis consists 
of three S functions, two P functions, and one D func- 
tion. Since the BFD pseudopotential removes a helium 
core for the first row atoms and a neon core for the sec- 
ond row atoms, the remaining elements lithium through 
argon have the same number of basis functions. In par- 
ticular, the 2z basis consists of two S functions, two P 
functions, and one D function, while the 3z basis consists 
of three S functions, three P functions, two D functions, 
and one F function. 

The bases consist of a combination of contracted and 
primitive functions. Since the BFD pseudopotentials are 
finite at the origin and have a Coulomb tail, the GS func- 
tions are the appropriate primitives. With the exception 
of the elements in Group \A of the periodic table (i.e. H, 
Li, and Na), the basis for each element includes a single 
S contraction and a single P contraction combined with 
an appropriate number of GS primitives. Only two con- 
tractions are employed to reduce the computational cost 
of using this basis in QC calculations. Since elements in 
Group \A of the periodic table have only one electron 



for the BFD pseudopotentials, a single S orbital is the 
ground state wavefunction, and this can be obtained ex- 
actly in HF. Thus, the basis for each element in Group 
\A includes a single S contraction, no P contractions, and 
an appropriate number of GS primitives. 



A. Contracted Functions 

A contracted basis function is a linear combination of 
Gaussian primitives: 
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where r, 6, (j) are the standard spherical coordinates, n 

is the principal quantum number, I is the azimuthal 
quantum number, m is the magnetic quantum number, 
Zp{0, (j)) is a real spherical harmonic, Cj is the i^^ expan- 
sion coefficient, and at is the i^^ Gaussian exponent. In 
practice, the restriction n = I + 1 applies. 

The exponents of the primitive functions that form 
the contracted basis functions are determined as follows. 
For each angular momentum for which a contraction is 
desired, an uncontracted basis consisting of nine even- 
tempered primitive Gaussians is generated. For each 
set of uncontracted Gaussians, the minimum exponent 
and even-tempering coefficient are varied to minimize 
the CCSD energy of the atom using a Python wrapper 
around GAMESS [ij. 

An assumption of weak coupling between the different 
angular momenta underlies the optimization procedure. 
Consequently, the uncontracted basis for each angular 
momentum is optimized separately. This optimization 
is performed by calculating the CCSD energy on an ini- 
tially coarse grid composed of different minimum expo- 
nents and even-tempering coefficients. Once regions of 
low CCSD energy are identified, a finer grid is used to ob- 
tain the final minimum exponent and even-tempering co- 
efficient. In addition to the assumption of weak coupling, 
two other properties of the problem make this global op- 
timization possible with modest computer resources; low 
dimensionality of search space and efficiency of atomic 
CCSD calculations. 

Next, an atomic MCSCF calculation in a complete ac- 
tive space (CAS) with the optimized uncontracted basis 
is performed in GAMESS. For these calculations, all elec- 
trons not removed by the pseudopotential are allowed to 
excite. For helium, the active space consists of the or- 
bitals from the n = 1 and n = 2 shells. For beryllium 
through neon, the active space includes the orbitals from 
the n — 2 and n = 3 shells. For magnesium through 
argon, the active space is composed of the orbitals from 
the n — 3 and n — A shells, with the exception of the 4D 
and 4F orbitals. A subset of the natural orbitals from the 
MCSCF calculations are used as the contracted functions 
of our basis. 



All atomic calculations are performed in D2h symme- 
try since GAMESS does not permit imposition of full 
rotational symmetry. Hence, different components of the 
same atomic subshell are not necessarily equivalent. Ad- 
ditionally, mixing may occur among orbitals of different 
angular momenta. For instance, there is mixing of S 



orbitals with both 03^2. 



and Y)x2_y2 orbitals. This 



anisotropy can be removed by averaging the different 
components of a particular subshell and zeroing out the 
off-diagonal blocks of the one-particle density matrix [7| . 
A simpler approach taken in this work is found to pro- 
duce results of similar quality. For each angular mo- 
mentum for which a contraction is desired, the NO with 
that angular momentum which has the largest occupa- 
tion number is chosen. Additionally, NO elements which 
do not correspond to the dominant character of the or- 
bital are zeroed out. For instance, an NO with large co- 
efhcients on the S basis functions and small coefficients 
on the D basis functions is considered to be dominated 
by S character, so the D coefficients are zeroed out. Fi- 
nally, the NOs are normalized. The NOs selected in this 
procedure generate the contracted functions for the basis 
set. The expansions of the contractions are given in the 
supplementary material [12l |. 



B. Gauss-Slater Primitives 

GS functions Q are defined as 

V^iimir, 0, 4>) = N< r"-ie-^ Z^e, cj,), (2) 

where C is the GS exponent and N^ is the normaliza- 
tion factor. The restriction n > I + 1 is imposed for GS 
functions. For r <C 1, the GS behaves like a Gaussian: 

'PiinAr, 0, 4>) = N^ r'-ie-^^'-)' Zr{9, 0), (3) 

and for r ^ 1, the GS behaves like a Slater: 

'fiiin.ir, 0, <t>) = Ni r"-ie-C- Z^C^, <t>)- (4) 

Consequently, GS functions introduce no cusp at the ori- 
gin and can reproduce correct long-range asymptotic be- 
havior of the orbitals. 

Unlike Gaussians and Slaters, normalization of GSs has 
no closed form expression. Nevertheless, normalizing an 
arbitrary GS is trivial with the following scaling relation 
between N^ and N}^: 






(5) 



Values for N^ are given in the supplementary material 

El. 

Since GSs are not analytically integrable, the radial 
part must be expanded in Gaussians for use in QC pro- 
grams that evaluate matrix elements analytically. The 



expansion is 


i 




^ii.ni'.H.'f)-- 





.c 



c 



(6) 
where c\ is the z*^ expansion coefficient and al is the 
i"^ Gaussian exponent. Notice that the expansion per- 
mits the case for which n ^ Z -f 1 for the GS function. 
Additionally, the following scaling relations hold for the 
expansion coefficients and Gaussian exponents: 
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Once the Gaussian expansions are found for unit expo- 
nents, expansions of arbitrary GSs follow immediately 
from the scaling relations. For QC calculations in this 
paper, GSs are expanded in six Gaussians. However, if 
the purpose of the initial QC calculation is to generate 
crude starting orbitals for QMC calculations in which 
orbital optimization is performed, it is only necessary to 
expand GS primitives in a single Gaussian. In this case, 
the cost of QC calculations is the same for Gaussian and 
GS primitives. The expansions of GS functions with unit 
exponent in both one and six Gaussians are given in the 
supplementary material [l^l ■ 

As mentioned above, the restriction n > I + 1 is im- 
posed for GS functions, instead of the more familiar 
n = I + 1 restriction imposed for Gaussian primitives. 
This motivates construction of two types of bases. In the 
first, ANO-GS, the restriction n = l+l is enforced. In the 
second, ANO-GSn, for each I there can be at most a sin- 
gle GS primitive with a particular n. For each additional 
primitive with a particular Z, n must be incremented. 

For example, consider lithium. The 2z ANO-GS basis 
has one S contraction, one GS-IS function, two GS-2P 
functions, and one GS-3D function. On the other hand, 
the 2z ANO-GSn basis has one S contraction, one GS-IS 
function, one GS-2P function, one GS-3P function, and 
one GS-3D function. 

A caveat to the above definition of the ANO-GSn ba- 
sis is that GS-2S functions are not permitted since a 
single GS-2S function will introduce an undesired cusp 
in the wavefunction. Additionally, the 2z ANO-GS and 
ANO-GSn basis sets are identical for all elements except 
lithium and sodium. When the 2z ANO-GS and ANO- 
GSn basis sets are identical, the basis sets are referred 
to as a 2z ANO-GS/GSn basis. For both lithium and 
sodium, the basis sets differ because these systems have 
no P contractions and instead have a second P primi- 
tive for the 2z basis. This primitive is a GS-2P for the 
ANO-GS basis and a GS-3P for the ANO-GSn basis. Ad- 
ditionally, weak coupling between functions of different 
angular momentum causes the GS-IS and GS-3D func- 
tions in the ANO-GS bases for lithium and sodium to 
differ from their counterparts in the ANO-GSn bases. 
However, the optimal exponents differ by less than 0.01. 



Optimal exponent selection for the GS primitives is dis- 
cussed now. Instead of optimizing exponents for the atom 
as was done to generate the contractions, optimization 
of the GS exponents is performed for the homonuclear 
diatomic molecule at experimental bond length jl4| - |22| . 
This advantageously produces a balanced basis set. 

Weak coupling between GS functions of different angu- 
lar momenta is assumed, so the initial optimization for 
each angular momentum is performed separately. This 
assumption is validated in Figure [1] which contains plots 
of the CCSD energy for Si2 while varying individual GS 
exponents in the 2z ANO-GS/GSn basis. Both the curve 
shape and exponent value which minimizes the energy 
vary little with fixed exponent value, signifying weak cou- 
pling between GS functions of different angular momen- 
tum. 
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FIG. 1. Change in Si2 CCSD energy for 2z ANO-GS/GSn 
basis shows weak coupling between GS functions of different 
angular momenta. TOP: Energy versus GS-IS exponent for 
three values of the GS-2P exponent with the GS-3D exponent 
fixed at its optimal value. Bottom: Energy versus GS-3D 
exponent for three values of the GS-IS exponent with the 
GS-2P exponent fixed at its optimal value. 



The optimization is performed at the CCSD level of 



theory using a Python wrapper around GAMESS. For 
each angular momentum, an energy landscape is defined 
by a grid of primitive exponents ranging from 0.1 to 6.0 
with 0.1 spacing. Thorough investigation has revealed 
that exponents larger than 6.0 are not optimal for the 
systems considered. Low lying minima of this energy 
landscape are then handled with increasingly finer grids 
until energy changes are less than 0.01 niH. During this 
investigation of local minima, all angular momenta are 
handled simultaneously to account for any coupling ef- 
fects. Results of this optimization are shown in Figure 
m Optimal exponents for ANO-GS and ANO-GSn bases 
exhibit a linear trend across each row of the periodic ta- 
ble. For nearly degenerate minima, the exponent follow- 
ing the trend in the figure is chosen as optimal, resulting 
in energy increase no greater than several 0.1 mH. The 
optimal GS exponents are given in the supplementary 
material [12j. 

In some cases, the optimal exponents for primitives 
with the same n and I are very close. This can lead to 
large equal and opposite coefficients on these basis func- 
tions when constructing molecular orbitals. Numerical 
problems could result, providing further motivation for 
the ANO-GSn basis, in which each pair of n and I is 
unique. However, all of our tests with the ANO-GS basis 
have had no numerical problems. 

Finally, the optimal primitive exponents are found to 
depend weakly on the electronic structure method em- 
ployed in the optimization, as demonstrated in Figure [3] 
for Si2 with the 2z ANO-GS/GSn basis. The globally 
minimizing exponents are nearly equal in different meth- 
ods. This exponent transferability to different levels of 
theory is extremely attractive for a basis set. 



III. RESULTS 

Section [II] demonstrates that the ANO-GS and ANO- 
GSn bases exhibit desirable properties. However, it re- 
mains to be shown that these basis sets produce accu- 
rate results. Fortunately, the basis set accompanying the 
BFD pseudopotential serves as a metric for testing ANO- 
GS and ANO-GSn basis quaUty. The BFD basis for el- 
ements in Groups lA and 2A of the periodic table has 
recently been updated [23| , but the number of functions 
in the new basis is inconsistent with the correlation con- 
sistent polarized basis prescription [l| . Since comparison 
would be difficult, their published functions are consid- 
ered in this work. 

Figure [5] shows the CCSD total energy gain per elec- 
tron of the ANO-GS and ANO-GSn bases over the BFD 
bases 10] for atoms and homonuclear dimers of hydro- 
gen through argon. Energy gains per electron tend to in- 
crease across each row of the periodic table. Both ANO- 
GS and ANO-GSn bases yield energy gains for most 
molecules and atoms. The energy gains per electron are 
generally larger for molecules than for atoms, and larger 
for the ANO-GSn basis than for the ANO-GS basis. The 
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FIG. 2. Optimal exponents for ANO-GS and ANO-GSn bases 
exhibit a linear trend across each row of the periodic table. 
The 2z ANO-GS and ANO-GSn bases are identical for all 
elements except lithium and sodium. The GS-IS and GS- 
3D exponents for these elements each differ by less than 0.01 
between 2z ANO-GS and ANO-GSn bases, so 2z ANO-GS 
and ANO-GSn are shown together as 2z ANO-GS/GSn. Ex- 
ponents for GS functions of P angular momentum are not 
included for lithium and sodium since these elements have an 
extra primitive of P angular momentum. 



FIG. 3. Change in Si2 energy for 2z ANO-GS/GSn basis 
shows optimal exponents depend weakly on electronic struc- 
ture method (CCSD, HF, and B3LYP). Top: GS-IS exponent 
is varied with GS-2P and GS-3D exponents fixed at their op- 
timal values. Middle: GS-2P exponent is varied with GS-IS 
and GS-3D exponents fixed at their optimal values. The large 
increase in energy around an exponent of 1.0 occurs since the 
P primitive and P contraction become nearly linearly depen- 
dent. Bottom: GS-3D exponent is varied with GS-IS and GS- 
2P exponents fixed at their optimal values. For Middle and 
Bottom, HF and B3LYP energy scale is on the right y-axis. 
This difference in energy scale occurs since higher angular 
momentum functions are less important in these effectively 
single-determinant theories. 



energy gains for the 2z bases are generally larger than for 
the 3z bases, as expected, since the energy left to recover 
becomes smaller as the basis size increases. 
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FIG. 4. CCSD total energy gains per electron of ANO-GS 
and ANO-GSn relative to the corresponding BFD basis [10] 
for atoms and homonuclear dimers of hydrogen through ar- 
gon. Energy gains per electron tend to increase across each 
row of the periodic table. The 2z ANO-GS and ANO-GSn 
bases are identical for all elements except lithium and sodium. 
Differences between 2z ANO-GS and ANO-GSn results for 
these elements is ~ O.OlmH, so they are shown together as 2z 
ANO-GS/GSn. 



The ANO-GS and ANO-GSn bases also produce more 
accurate CCSD atomization energies than the BFD basis 
for the homonuclear dimers of hydrogen through argon. 
Figure [S] shows the fraction of experimental atomization 
energy recovered in CCSD for the homonuclear dimers 
which are not weakly bound. The 2z ANO-GS/ANO- 
GSn basis recovers more atomization energy than the 2z 
BFD basis for all dimers except those of Group lA ele- 
ments. Similarly, the 3z ANO-GSn basis recovers more 
atomization energy than the 3z BFD basis for the same 
systems, but the differences are small. The 3z ANO-GSn 
is on average slightly better than the 3z ANO-GS basis, 
the largest gains being for F2 and CI2. 

For Group lA elements, the BFD bases recover more 
atomization energy in CCSD than do their ANO-GS or 
ANO-GSn counterparts. This occurs due to inaccurate 
BFD energies for the atoms, as can be seen in Figure 
m However, as described above, we used the published 
BFD bases for these elements rather than the updated 
BFD bases 23] to maintain consistency. 

Finally, improvements of the ANO-GS and ANO-GSn 
bases extend to other systems and methods. Figure [S] 
shows the fraction of experimental atomization energy 
recovered for five systems in the G2 set [2y] with the 
BFD, ANO-GS, and ANO-GSn bases in three quantum 
chemistry methods. For CCSD, the ANO-GS and ANO- 
GSn bases outperform the BFD basis for all systems. For 
sulfur dioxide the improvement due to the ANO-GS and 



ANO-GSn bases is dramatic: the 2z ANO-GS/GSn result 
is nearly halfway between the 2z and 3z BFD results, and 
the 3z ANO-GS/GSn result is nearly halfway between the 
3z and 5z BFD results. ANO-GS and ANO-GSn benefits 
are more prominent in HF and B3LYP: for most systems, 
the 2z ANO-GS/GSn resuft is closer to the 3z BFD resuft 
than the 2z BFD result, and the 3z ANO-GS/GSn resuft 
is closer to the 5z BFD result than the 3z BFD result. 
Differences between results with the ANO-GS and ANO- 
GSn bases are small. 

Figure [7] shows the fraction of experimental atomiza- 
tion energy recovered using diffusion Monte Carlo (DMC) 
with the BFD, ANO-GS, and ANO-GSn bases. For each 
system, the DMC calculations are performed with both 
a single-configuration state function (single-CSF) refer- 
ence (DMC-ICSF) and full- valence complete active space 
reference (DMC-FVCAS). However, for each of the con- 
stituent atoms in these molecules, the FVCAS and single- 
CSF references are equivalent. All DMC calculations are 
performed with a 0.01 H"-'^ time step and trial wavefunc- 
tion obtained by optimizing Jastrow, orbital, and config- 
uration state function (CSF ) pa rameters (where applica- 
ble) via the linear method |28l43(j | in variational Monte 
Carlo. The DMC-ICSF and DMC-FVCAS calculations 
exhibit similar trends to the HF and B3LYP calculation 
for most systems: the 2z ANO-GS/GSn result is closer 
to the 3z BFD result than the 2z BFD result, and the 3z 
ANO-GS/GSn resuft is closer to the 5z BFD result than 
the 3z BFD result. Again, differences between results 
with the ANO-GS and ANO-GSn bases are smah. 

There are several important points that can be made 
by comparing the DMC calculations of Figure [7] to the 
CCSD calculations of Figure M Fhst, the DMC results 
for the atomization energies have a weaker dependence 
on basis size than the CCSD results. Second, for a given 
basis set, the most basic DMC calculations, DMC-ICSF, 
yield superior results compared to CCSD. In addition to 
yielding superior results, DMC-ICSF calculations have 
better computational cost scaling than CCSD calcula- 
tions. Under certain assumptions, the cost of DMC-ICSF 
calculations scales as 0{N^) l3lt, while the cost of CCSD 
calculations scales as 0{N^) [S^l, where A'' is the num- 
ber of electrons. However, it is important to note that 
the prefactor of the scaling is significantly smaller for the 
CCSD calculations. 

Finally, our results are not the first to show that DMC 
calculations can produce accurate atomization energies. 
In particular, DMC-ICSF calculations of the entire G2 
set have been performed for both pseudopotential and 
all-electron systems [33, |3J| and produced excellent re- 
sults. Additionally, there is good agreement between the 
pseudopotential and all-electron results with a mean ab- 
solute deviation of about 2.0 kcal/mol over the entire G2 
set [3J]. Although these previous results are very good, 
there is room for improvement, particularly for the open 
shell systems. A systematic study with DMC-FVCAS 
calculations is currently underway in our group, which 
should produce results to (near) chemical accuracy for 





1.2 






CI) 




d 


1 1 


LU 




D) 


1 


c 


"O 




c 


0.9 


DQ 




d. 


0.8 


X 




LU 


0.7 


O 




o 


0.6 


'.^^ 




!^ 


0.5 


i— 




II 






0.4 




FIG. 5. Fraction of experimental atomization energy recovered in CCSD with BFD, ANO-GS, and ANO-GSn bases for the 
homonuclear dimers of hydrogen through argon which are not weakly bound. The 2z ANO-GS and ANO-GSn bases are 
identical for all elements except lithium and sodium. Differences between 2z ANO-GS and ANO-GSn atomization energies for 
these elements is ~ O.OlmH, so they are shown tog ether as 2z ANO-GS/GSn. Calculated values are corrected for zero point 
energy [13, [Ij to compare with experiment p3. [ig. [TtI. [25l] . 



all systems in the G2 set. 



IV. CONCLUSION 

A simple yet general method for constructing basis sets 
for molecular electronic structure theory calculations has 
been presented. These basis sets consist of a combination 
of atomic natural orbitals from an MCSCF calculation 
with primitive functions optimized for the corresponding 
homonuclear dimer. The functional form of the primitive 
functions is chosen to have the correct asymptotics for the 
nuclear potential of the system. 

It was shown that optimal exponents of primitives with 
different angular momenta are weakly coupled. This en- 
ables efficient determination of optimal exponents. Ad- 
ditionally, it was demonstrated that the particular elec- 
tronic structure method employed in optimization has 
little effect on the optimal values of the primitive expo- 
nents. 

Two sets of 2z and 3z bases, ANO-GS and ANO-GSn, 
appropriate for the Burkatzki, Filippi, and Dolg non- 
divergent pseudopotentials were constructed for elements 
hydrogen through argon. Since these pseudopotentials 
do not diverge at nuclei and have a Coulomb tail, GS 
functions are the appropriate primitives. 

It was demonstrated that both ANO-GS and ANO- 
GSn basis sets offer significant gains over the Burkatzki, 
Filippi and Dolg basis sets for CCSD, HF, B3LYP 11], 
and QMC calculations. Improvements were observed in 
both total energies and atomization energies. The latter 
indicates that basis sets providing a balanced description 
of atoms and molecules were produced by using both the 
atom and the dimer in the optimization. On average, 



the ANO-GSn basis is slightly better than the ANO-GS 
basis, but either is a sound choice. 

In the future, these basis sets will be extended to 
include the transition metals, and, bases will be con- 
structed for all-electron calculations, for which Slater 
functions are the appropriate primitives. 
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FIG. 7. Fraction of experimental atomization energy recov- 
ered in difltusion Monte Carlo (DMC) for LiF, O2, P2, S2, and 
SO2 with the BFD, ANO-GS, and ANO-GSn bases. DMC 
calculations are performed with both a single-CSF reference 
(DMC-ICSF) and full- valence complete active space reference 
(DMC-FVCAS). The 2z ANO-GS and ANO-GSn bases yield 
different results only for LiF. The 5z BFD* calculations do 
not include the G or H functions from the 5z BFD basis. 
Calculated atomization energies are corrected for ze ro p oint 
energy [13, [111 to compare with experiment [ll. [l^. [TtI. [27J . 
The legend for this plot is identical to that of Figure |6l 
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FIG. 6. Fraction of experimental atomization energy recov- 
ered in HF, B3LYP, and CCSD for LiF, O2, P2, S2, and SO2 
with BFD, ANO-GS, and ANO-GSn bases. The 2z ANO- 
GS and ANO-GSn bases yield different results only for LiF. 
The 5z BFD* calculations do not include the G or H func- 
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